/*
 * mpole_misc.c - miscellaneous subroutines for multipole library
 *
 * w. elliott
 * w. t. rankin
 *
 * Copyright (c) 1995 Duke University
 * All Rights Reserved.
 *
 */

/*
 * include files
 */
#include "mpole_misc.h"

     

/****************************************************************
 *
 *  dumpY_C() - print out a multipole expansion to standard error.
 *
 */

void dumpY_C( Mtype Y, int p )
{
   int n, m;

   for (n = 0; n < p; n++) {
      for (m = 0; m <= n; m++) {
	 fprintf(stderr, "%.6e %.6e   ", Y[n][m].x, Y[n][m].y);
      } /* for m */
      fprintf(stderr, "\n");
   } /* for n */
   fprintf(stderr, "\n\n");

} /* dumpY_C */


/****************************************************************
 *
 *  dumpYF() - print out an FFT multipole expansion to standard error.
 *
 */

void dumpYF( Real *Ystart, int p )

{
   int n, m, pf;

   pf = 1 << ((int) (log((double) (2 * p - 1)) / log(2.0)));

   for (n = 0; n < 2 * p; n++) {
      for (m = 0; m < pf; m++) {
	 fprintf(stderr, "%e\t%e\n", Ystart[2*n*pf+2*m], Ystart[2*n*pf+2*m+1]);
      } /* for m */
      fprintf(stderr, "\n");
   } /* for n */

} /* dumpYF */


/****************************************************************
 *
 *  MathdumpY_C() - dumps a multipole expansion to a file
 *    specified in s[].
 *
 */

void MathdumpY_C( Mtype Y, int p, char s[], int cell )
{
   int             n, m;
   FILE           *fpo;

   fpo = fopen(s, "a");

   fprintf(fpo, "mpole%d = [\n",cell);

   for (n = 0; n < p; n++) {
      for (m = 0; m < p; m++) {
	 if (m <= n) {
	    fprintf(fpo, "%g", Y[n][m].x);
	    if (Y[n][m].y > 0.0)
	       fprintf(fpo, "+i*%g ", Y[n][m].y);
	    else if ( Y[n][m].y < 0.0 )
	       fprintf(fpo, "-i*%g ", -Y[n][m].y);

	 } /* if valid array value */ 
	 else
	    fprintf(fpo, "0");
	 if (m != p - 1)
	    fprintf(fpo, " ");
      } /* for m */
      fprintf(fpo, ";");
      if (n < p - 1)
	 fprintf(fpo, " ");
      fprintf(fpo, "\n");
   } /* for n */
   fprintf(fpo, "]\n\n");

   fclose(fpo);

} /* MathdumpY */

/****************************************************************
 *
 *  MDumpRaw_C() - dumps a multipole expansion to a file
 *    specified in *s.
 *
 */

void MDumpRaw_C( Mtype Y, int p, char *s )
{
   int             n, m;
   FILE           *fpo;

   fpo = fopen(s, "w");

   for (n = 0; n < p; n++) {
      for (m = 0; m <=n; m++) {
	 fprintf(fpo, "%20.16g %20.16g\n", Y[n][m].x, Y[n][m].y);
      }
   }
   fclose(fpo);

} /* MathdumpY */

/****************************************************************
*
*  Cart2Sph(v,s) - converts cartesian to spherical vectors
*
*/
void Cart2Sph( vektor v, vektor *s )
{
   Real x, y, z;
   Real rho, alpha, beta;

   x = v[0]; y = v[1]; z = v[2];

   rho = sqrt(x*x + y*y + z*z);

   if (rho < fabs(z)) rho = fabs(z);

   if (rho == 0.0) alpha = 0.0;
   else alpha = acos(z/rho);

   if ((x == 0.0) && (y == 0.0)) beta = 0.0;
   else beta = atan2(y, x);

   (*s)[0] = rho; (*s)[1] = alpha; (*s)[2] = beta;

} /* Cart2Sph */
